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We review recent statistical mechanical approaches to multimode laser theory. The theory has 
proved very effective to describe standard lasers. We refer of the mean field theory for passive mode 
locking and developments based on Monte Carlo simulations and cavity method to study the role of 
the frequency matching condition. The status for a complete theory of multimode lasing in open and 
disordered cavities is discussed and the derivation of the general statistical models in this framework 
is presented. When light is propagating in a disordered medium, the system can be analyzed via the 
replica method. For high degrees of disorder and nonlinearity, a glassy behavior is expected at the 
lasing threshold, providing a suggestive link between glasses and photonics. We describe in details 
the results for the general Hamiltonian model in mean field approximation and mention an available 
test for replica symmetry breaking from intensity spectra measurements. Finally, we summary some 
perspectives still opened for such approaches. 


The idea that the lasing threshold can be understood as a proper thermodynamic transition goes back since the 
early development of laser theory and non-linear optics in the 1970s, in particular in connection with modulation 
instability (see, e.g.^ and the review^). A proper formalization of this idea was achieved in the ’00s, with the works of 
B. Fischer et ali^“— , where the statistical properties of laser light in homogeneous cavities are investigated taking into 
account nonlinear effects, like gain saturation and intensity dependent refractive index. Mapping the laser dynamics 
in an ordered Hamiltonian problem, with the nonlinearities providing the interaction between the electromagnetic 
modes, the authors of Reffi show new aspects and reveal, in particular, that a critical behavior occurs at the laser 
mode locking transition. In section U of this paper we review the main results and the most recent developments^^— 
of this use of statistical mechanics of ordered systems for the study of light propagation in nonlinear materials in 
homogeneous cavities. 

In the last years the role of local inhomogeneities in matter and its disordered nature has proved to be at the origin 
of a variety of novel interesting phenomena. Specifically, light amplification in random media and random lasers 
(RL) have attracted much attentionii^— . This is an intriguing topic that crosses different fields as light localization, 
non-linear physics and quantum optics and has relevant fundamental aspects and as well practical perspectives like, 
e.g., biomedical diagnostics^!, chip-based spectrometers!!"— and cryptography!!. 

From a historical point of view, the presence of a divergence in the intensity for light diffusion with gain above 
a critical volume was already discussed by Letokhov in the 1960s. However, lasing occurs in general as a result of 
two basic ingredients: optical amplification and feedback. Amplified spontaneous emission occurs even without optical 
cavity, the spectrum is completely determined by the gain curve of the active material. The material is called a 
RL only when the multiple-scattering process (feedback) plays a key role in determining the lasing process!!. The 
presence of feedback is associated to the existence of well-defined cavity modes with a long-life and characterized by a 
dehnite spatial pattern of the electromagnetic field, which sets up inside the lasing structure in the stationary regime. 

A basic theory for RLs may, in principle, be developed along the lines of standard semiclassical lasing theory, but 
it has to incorporate the specific features of random lasers: a much larger role of radiative leakage and an irregular 
spatial structure of the modes. While this theory is far from being complete, a review of some recent developments, 
including the essential elements for a statistical approach, is reported in section [TTl In this physical situation, the 
methods of statistical mechanics have proved to be very powerful to treat the strong interplay between disorder and 
nonlinearity!!^— . Within this approach, it is possible to draw suggestive analogies between propagation of light in 
nonlinear disordered media and the rich behavior observed in glasses. We report some of the most important results 
from this approach in section uni 


I. THE CASE FOR STANDARD LASERS 

The statistical mechanics paradigm called Statistical Light-mode Dynamics (SLD) has been introduced in the early 
works!i! to study multimode laser physics, for which the number of modes is generally large enough (10^ — 10® in 
long lasers), and such that the presence of nonlinearity makes the problem nontrivial. One of the main benefits 
of this approach is the recognition of the role of the noise, due, e.g., to the unavoidable spontaneous emission. In 


2 


literature, see, e.g., RefS i^^i^° , classical and quantum noise is often considered as a small perturbation of the noiseless 
pulse of the laser master equation. In the SLD approach, instead, the noise takes the role of the temperature in the 
thermodynamic system and can be treated nonperturbatively. As a particular result, it is shown that the perturbative 
approach is in general not legitimate: accumulation of noise in the whole cavity can generate a continuous background 
that carries a significant portion of the total optical power and competes with the pulsed. 

The entropy associated with the noise becomes an essential ingredient in the study of mode locking, the regime 
under which a laser generates ultra-short pulses^. The pulse formation in lasers is based on the interaction between 
axial modes. Such an interaction can be provided either by making the system time dependent {active mode locking) 
or by a suitable nonlinearity in the dynamics of the system {passive mode locking)^^ Within the SLD the difference 
between active and passive mode locking is clear and it is embedded in the range of the interaction between modes in 
the two cases. The active case corresponds to a one-dimensional short-range-interacting model^: in this case a phase 
transition to an ordered state occurs in principle only at zero noise. The fragility of active mode locking becomes, 
hence, an exemplification of the well-known lack of global ordering of the one-dimensional spherical spin model^: any 
weak noise can break a bond between two interacting modes, thus eliminating global mode ordering. In contrast, 
the long-range interaction in passive mode locking, due to the four-wave mixing in the saturable absorber, imposes a 
global order below a certain noise level, resulting in the threshold behavior: once the interaction is strong enough to 
overcome noise it induces long-range correlations and a first order phase transition between disordered (continuous 
wave) and ordered (mode locked) thermodynamic phases^‘^‘^i2&. Thus, the individual action of the noise stabilizes the 
continuous wave regime, showing a noise induced phase transition^. 

Several other theoretical and experimental features can be studied in SLD. Among the others, we mention the 
addition of an external driving field to a passive mode locked system, analog to the external magnetic field in magnets 
or to the pressure in gas-liquid-solid systems, representing the injection in the laser cavity of pulses from an external 
source, which in the simple case matches the repetition rate of the laser—. In this situation, beyond a threshold 
injection level, the transition becomes continuous rather than first order. The two phase transition lines meet at a 
tricritical point around which tricritical behavior is observed. 

For a passive mode locking laser, the SLD description can be obtained considering the standard master equation^ 

^ = {Gi + iDi) ai{t) + {T - iA) ^ ak,{t) ■ al^{t) ■ ak,,{t) + Fi{t), (1) 

fcl—^‘2+^3 —^ —0 

where the slow amplitudes a/ at the frequency uji = uiq +1 Aui, as imposed by the Fabry-Perot resonator, are defined 
via the expansion of the electric field 


E{z, t) = Y, -b c.c.. (2) 

i 

In Eq. ((T|) the real parameter Gi represents the difference between the gain and loss of the mode I in a complete 
round-trip through the cavity, Di is the group velocity dispersion of the wave packet, F is the nonlinear self-amplitude 
modulation coefficient associated to a saturable absorber and A is the self-phase modulation coefficient (responsible 
of the Kerr lens effect). The noise Fi{t) is generally assumed Gaussian, white and uncorrelated: 


{n,{ti)FkAt2)) = 2T5k,kJ{ti-t2), {FkAti)Fk,{t2)) = 0, (3) 

being T the spectral power of the noise, associated, e.g., to spontaneous emission. 

To obtain an effective statistical approach the gain saturation must be explicitly considered. Indeed, when Gi and 
F are nonzero, the total optical intensity E = ^ constant of motion for the evolution of Eq. ([1]). In 

this case the laser is stable because the gain decreases as the optical intensity increases^^. For standard lasers this 
is usually modeled assuming that the gain in the master equation Eq. o is given by G; = Go/(I +£/Esat), being 
Esat the saturation power of the amplifier. To study the equilibrium properties of the model, in SLD it is assumed a 
simpler model: at any instant the gain is supposed to assume the value that keeps £ exactly a constant of motion^. 
This corresponds to assuming the spherical constraint £ = Eq during the evolution of the slow amplitudes ai. The 
relation between the thermodynamics in the fixed-power ensemble and a variable-power ensemble is similar to the one 
between the canonical and grand canonical ensembles in statistical mechanics^. 

In the limit of small dispersion Di <C Gi and A <C F, the evolution becomes Hamiltonian with 


% = - Y. 

k=\,...N 


F 
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Ofci ttk^ ak^ ak^ , 

fe 1 —/l‘2 + ^3 —/C4 =0 


( 4 ) 


where, in the fixed-power ensemble, Gi is constant and the condition £ = Eq = eN is preserved. This model can be, 
then, directly investigated with standard methods of statistical mechanics for a system at the effective temperature 
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Tph = T being V the pumping rate of the source. In particular, the mean field theory (MFT) of the model 
is exact when the Frequency Matching Condition (FMC) fci — fc 2 + ^3 ~ ^4 = 0 on the nonlinear term is neglected. In 
this case the above mentioned first order transition is obtained between a continuous and a pulsed regime. 

The role of the FMC in mode locked laser systems has been the object of recent works within the statistical 
approach2r2£. It is observed that, even in a laser system in the presence of FMC, a first order transition separates 
a continuous wave regime at low pumping and a pulsed regime at high pumping, the critical pumping value being 
compatible in the thermodynamic limit with the MFT one. The difference between the MFT and FMC solutions 
is in the nature of the pulsed regime. In MFT, the modes are all equivalent and they are trivially phase-locked 
with the same phase (pj = (po for j = 1,...N at high pumping. Alternatively, when the FMC is included the 
pulsed regime is not trivially phase locked and the laser mo de ph ases result to behave as pj ~ pQ + p' uj, with a 
nontrivial, frequency independent, slope p' (see insets in Fig. |l(a)[ ). The slope p' changes in the time evolution with 
a distribution determined by the interaction network among the modes and with a lifetime that increases with the 
pumping rate^d^. Consequently, all the two mode correlators are zero in the high pumping regime in FMC system&i^. 
The electromagnetic pulse associated with such configurations, cf. Eq. 0, is, accordingly, such that the phase delay 


between the carrier and the envelope of the signal changes at each shot, cf. Fig. 1(a) 


The FMC turns out to play an essential role also for the intensity spectra: for low pumping the spectrum simply 
follows the gain curve Gk of Eq. o, but, as the pumping exceeds the mode locking threshold, the spectrum is 
determined mainly by the interaction network. In this case, then, it is observed to become narrow around the central 
frequencies of the bandwidth, cf. Fig. |l(b)[ as a direct consequence of the inhomogeneous connectivity induced by the 
FMC. This furnishes a simple theoretical mechanism to explain the gain narrowing at the mode locking transition. 

The transition from a continuous regime to a pulsed regime is observed also in models with fixed intensities^d^: the 
fluctuations of the amplitudes is neglected with respect to the one of the phases and the resulting Hamiltonian ([3]), 
obtained rescaling the interaction coefficients with the averaged amplitudes, can be mapped to the problem of non- 
linearly interacting XY spins, a = {cos p, ship). In Refi^, the equilibrium configurations of the phases are obtained 
in the thermodynamic limit mapping the system on random regular graphs: when the number of modes is very high, 
the details of the graph become unimportant and the solutions are found averaging over many realizations of the 
network. The equilibrium configurations of the system are then studied considering an extensive number of non- 
linearly interacting quadruplets. This limit is interesting for RLs in which the spatial distribution of the modes could 
be highly localized. Instead to employ the Monte Carlo method, one can then study these diluted network systems 
using the cavity method, widely used for optimization problems with discrete variable nodes. The cavity method allows 
to determine the equilibrium distributions of the variable nodes starting from marginal probability distributions: to 
each link connecting a variable node and a function node are associated two probability distributions, also called 
messages, representing the marginal probability distribution of one variable when the other is not there; e.g., r/i^a{P) 
indicates the message sent from variable node i to function node a. The idea is that, since the graph is highly diluted, 
when a link is cut the other variable nodes participating in that function node become uncorrelated; the problem can 
then be solved self-consistently: 


where with h G di\a we indicate all the neighbors of i except a. Irregularity in the graph are considered evaluating 
the distribution of the rjs. Non-trivial equilibrium configurations are obtained for low values of the effective bath 
temperature: together with the ferromagnetic phase, in which all the modes oscillates synchronously and that breaks 
the 0(2) symmetry, a phase wave is observed. The electromagnetic field is also reconstructed from the (pj) and 
the different regimes are analyzed also respect to the starting frequency distributions. An electromagnetic pulse is 
obtained for a frequency comb distribution. 


II. MULTIMODE LASER THEORY FOR OPEN AND DISORDERED LASERS 

After the pioneer works of Letokhov^, there was a renewed interest starting from the 90’s in light amplification in 
disordered systems both for fundamental and applied research. The main features that distinguish RLs with respect to 
conventional lasers are random disorder in the refractive index and absence of a light trapping cavity: the confinement 
is realized by the multiple-scattering caused by the disorder itself. In the first experiments it was unclear whether 
even this cavityless systems could have well-defined electromagnetic modes showing narrow peaks in the spectrum. 
Indeed, light intensity was described by a diffusion equation, with phases and hence interference not playing any role. 
In these cases, the threshold for lasing is obtained when the amplification length becomes smaller than the path length 
a photon travels before escaping from the gain medium. Instead of high-Q resonances characterized by discrete peaks 
in the spectrum, as in more conventional lasers, there appear low-Q resonances spectrally overlapping: the spectrum 
is continuous, it does not contain any discrete components at selected resonance frequencies and it is centered at 
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(a) Signal of Eq. (top) and corresponding 
configurations plotted as Phase vs Frequency (bottom) 
for a laser system Eq. 0 with 500 modes just above 
the critical pumping {V = 1.692, Vc = 1.565(8)). 
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(b) Intensity spectra for three pumping for a laser 
system Eq. Q with 150 modes {Tc = 1.60(2) here). 
Below the threshold the spectrum is shown to follow 
the gain curve G{cok) = Gk (black solid line). 



the atomic transition, which is the only remaining resonant element of the system. The feedback in these systems is 
indicated as intensity or energy feedback, it is incoherent and non-resonant (for a review on non-resonant feedback 
random lasers see Ref.—). 

In later experiments^, focusing higher intense pumps in smaller systems, it was possible to actually observe discrete 
lasing peaks in the intensity spectrum. These peaks were not necessarily centered around the atomic gain line and 
exhibited a Poisson photon count distribution, typical of coherent lasing light. Moreover, the position of the peaks 
was observed to change when a different part of the sample was illuminated. It was suggested that this transition in 
the spectrum was due to a transition from a nonresonant to a resonant feedback: even though RLs are characterized 
by the absence of a well-defined cavity yielding strong leakages, electromagnetic modes have enough time to form and 
lasing involves coherent phase-sensitive feedback. 

After these new results, it has been accepted that RLs are multimode systems and research in this field has been 
focusing on the development of a multimode laser theory, which, starting from the lines of a standard semiclassical 
multimode laser theory^, has to incorporate specific features of these systems: strong radiative leakage, modes with 
a broad distributions of lifetimes and irregular spatial dependence of the refractive index; for three dimensional (3D) 
systems this implies also the absence of well-defined optical cavity axes. 

In the rest of this section we will briefly reference several alternative ways to introduce a multimode descriptions 
of strongly open systems. Starting from a quantum treatment, we will then analyze in more details the approach of 
Hackenbroich and co-worker a^'^d^ obtaining, under some hypothesis, a Hamiltonian of the form of Eq. O- 


Cold cavity modes in highly open systems 

The problem of a modal description for open systems is of interest not only for RL physics and there have been 
introduced many descriptions for this purpose. The difficulties of the problem increase due to the fact that leakages 
make the system non-Hermitian and the usual technique of developing a basis of orthogonal eigenvectors with real 
eigenvalues does not apply. Concerning laser physics. Fox and Li^ firstly studied the effect of losses for the electro¬ 
magnetic modes in a cavity. In particular, they showed how the low-order modes, i.e. the modes that have lowest 
threshold, could be determined by a numerical iterative method: they assumed an arbitrary field distribution at the 
first mirror and proceed to compute its distribution after one passage in the cavity; this distribution is then used to 
compute the field at the second mirror and so on. The modes are then determined by the steady state distributions 
solving the self-consistent equations: 


/ A'(r_L,r'_L)^„(r'_L)dr'_L = A„V'ri(i*j.) (5) 

J A 

where A is the mirror surface, K is the kernel determined through the Huygens-Fresnel formulation and T indicate the 
direction transverse to the optical axis; A„ is a complex constant whose module and phase determine the attenuation 
and phase shift per transit for the mode. The application of Fox-Li modes for random lasing physics is rather limited 
for several reasons. First, in the 3D case, RLs lack of a well-defined optical axis and, more in general, the convergence 
of the self-consistence equations, Eq. becomes slow and less effective when there are a lot of modes, which is the 
typical case for RLs. 
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In general, the approach adopted in laser physics to predict the lasing modes is to start from the cold cavity 
modes, which are solutions of the Maxwell equations with proper boundary conditions. In order to describe leakages, 
there have been introduced the so-called quasibound states^: solutions of the Maxwell equations without gain and 
with boundary conditions containing only outgoing waves at infinity. The Fox-Li method previously described is a 
technique to find the quasibound states of the system. The quasibound states are, though, only defined inside the 
cavity: they are associated to complex eigenvectors, see Eq. and, if extended outside the cavity, diverge at infinity, 
since they satisfy Maxwell equations with only outgoing wave. For this reason, quasibound modes cannot be used for 
a quantization procedure in the entire space. 

The standard procedure to go beyond the cold cavity modes is to use the semiclassical laser theory, including the 
gain medium described as a two levels atom. The atoms are pumped to the upper level and then can decay to the 
lower level emitting a photon in the mode field; other atomic levels are in general introduced through a finite lifetime 
of the two laser levels. The pumped medium is coupled to the field satisfying the Maxwell equations. The equations 
obtained are the Maxwell-Block equations (see, e.g., Ref4^ for the case of a cavity with a well-defined optical axis). 
One usually finds, for high-Q cavities, that the mode with highest Q lases first with a frequency shifted toward the 
atomic transition. The spatial distribution of the mode is not supposed to be altered much by the gain medium, i.e., 
it is basically the same distribution of the cold cavity mode. The difficulties of the Maxwell-Block equations comes 
from the presence of multimodes: in this case nonlinearities can develop in the system and the equations cannot be 
solved exactly. Furthermore, even disregarding nonlinearities, a legitimated hypothesis only below the threshold, the 
interactions are nontrivial in presence of many modes since they will depend on the (complex) modes spatial overlap. 

In respect to high-Q cavities, for the case of highly open systems with a nonuniform refractive index, the lasing modes 
can be quite different from the cold cavity modes. One approach that has been used in this case is to determine the 
spatial distributions of the lasing modes as well as their frequencies solving Maxwell-Block equations self-consistently. 
In Ref4i the authors apply this approach disregarding nonlinearities caused by possible beat frequencies, i.e. neglecting 
population pulsation. To analytically obtain the solution, they use the Green function written in terms of a new set of 
eigenfunctions, the constant flux states. Constant flux states are similar to the quasibound states, having a complex 
wave vector inside the cavity, but differ from the latter because they have a real wave vector outside the cavity 
not carrying a diverging flux at infinity (whence the name “constant flux”). Constant flux states solve different 
differential equations inside and outside the cavity; continuity determines boundary conditions. Considering the case 
of a one-dimensional slab, the set of the constant flux states {<j>m{x,uj)} is defined to satisfy 

-V^^m(a;,a;) = n^(a;,w) {km + kaf 4>m{x,uj) (6) 

inside the cavity, while 

-^'^<l>m{x,u:) = {k + kaf 4>m{x,u}) (7) 

outside the cavity with outgoing wave boundary conditions. Here n(a;, w) is the refraction index of the material and 
ka = Ea {hc)~^, being Ea the atomic transition energy. Constant flux states are not orthogonal or complete. Due 
to the non-Hermitian nature of the problem, to obtain a complete set, it is also necessary to introduce the functions 
{<j>m}j which are the adjoint of the constant flux states. These correspond to the solution for free propagation outside 
the cavity with incoming waves only. Constant flux states with their adjoint functions form a complete set that can 
be used for a quantization in the whole space: the amplitudes become operators, and two set of operators are used 
inside and outside the cavity. However, constant flux states are not very convenient as a basis for a quantization 
procedure since the internal and the external operators do not commute; non-commuting relations introduce coupling 
between the inside and the outside operators. For this reason, in the rest of this section, we will introduce the 
system-and-bath-approach as a starting point for a quantization procedure. 


The system-and-bath approach 

The idea of the quantization technique based on the system-and-bath-approach is to separate the space in two 
subsystems introducing separate operators and eigenmodes. One subsystem is the open cavity, the resonator, con¬ 
taining the gain medium; the open surrounding space constitutes the bath. The bath and the resonator interact at 
the boundary. This technique is adequate to study the field degrees of freedom in the resonator: at the end of the 
quantization procedure the degrees of freedom of the bath can be averaged out, obtaining then information about the 
field inside the resonator alone. 

The system-and-bath model is in general a standard description of losses in quantum mechanical system o^^d^ and 
it has been used for more than 40-years for quantum optical systems. Some well-known, textbook examples, are 
the derivations of the finite lifetimes of the excited levels of an atom and the damping of a “single” oscillator mode. 
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Less understood is instead the damping in the case of multimode fields. In this case, the system is modulated by a 
discrete set of independent harmonic oscillators, obtained from the quantization of the normal modes of the system. 
A continuum set of harmonic oscillators models the outside reservoir. Coupling between the two gives rise to damping. 
However, in this description, coupling constants only enter as phenomenological parameters and it has been argued 
that it can be applied only to systems with low coupling with the outside reservoir, as more conventional lasers. 
Hackenbroich and coworkers have shown in a series of papers^^— how to rigorously derive the system-and-bath 
Hamiltonian starting from Maxwell equations employing the Feshbach projector formalism^, a well-known technique 
in condensed matter and nuclear physics, hence providing a microscopic description of the theory. In Ref4i, they also 
derived explicitly the expression for the coupling constants for some simple geometries. We will briefly sketch their 
results. 

Let us consider the case of uniform three-dimension linear dielectric medium, e(r). The medium is surrounded by 
free space, e(r) = 1. In the absence of source, the eigenmodes fm(<j-’,r) of the whole system satisfy: 

V X [V X fm(a;,r)] - = 0 (8) 

the fm(w,r) are defined in the whole space, w is a continuous parameter while the discrete index m indicates the 
boundary conditions depending on the problem at hand. For example, for a dielectric coupled to free space m indicates 
the possible polarizations. In order to transform Eq. (I8|) into an eigenvalues problem, it is convenient to introduce 
the set of functions, r) = A/e(r)fm(w, r); the r) are the solution of the equation 


L(j){uj,r) = 


\/eW 


V X 


V X 


\/eW 




( 9 ) 


The functions (pm{ui,r) form a complete set in the subspace of functions defined by the condition V • 
e(r)(/>m(w, r) = 0. The vector potential A can be expressed as a superposition of the (^m(w,r) functions with 
proper coefficients 


A(r,t) = cy ; 


( 10 ) 


same for the canonical momentum, n(r) = e(r)A(r)/c^, with coefficients Quantization is now obtained 

promoting q^uj) and p{u!) as operators with proper commutation relations^. 

The quantization of the electromagnetic modes starting from the eigenfunctions of the Helmholtz equations, Eq. 
(|8|). is known as the mode-of-the-universe approach. The resulting Hamiltonian is a sum of harmonic oscillators. 
However, as discussed before, we would like to have explicit information on the modes inside the resonators. In the 
system-and-bath approach, we consider then two subsystems: one containing the resonator the other the open space 
outside. This can be implemented defining the two Feshbach projection operators^ 


Q = [ dr|r)(r|, ^ = / <^i‘|r)(r| (11) 

Jr£V Jr^V 

representing the projection into the resonator and in the outside space respectively; V indicates the space occupied 
by the resonator. The operators Q and V satisfy: Q = Q^Q = Q^, same for V, and are complete and orthogonal 
Q + V =1, QV = VQ = 0. The functions can now be decomposed through V and Q 


<)'(r) = Q<)'(r) + 'P(f){r) = /i(r) -|- v{r) 

and Eq. m can be expressed, then, as 

fL qq / p{r)'\ ^ 

\L-pq L-pvJ J \cJ [i^(r)J 


( 12 ) 


(13) 


with Lqq = QLQ, Lqp = QLV and so on. The operators Lqq and L-p-p are self-adjoint operators in the Hilbert 
space of the resonator and outside channel functions, respectively. Let us indicate the set of functions satisfying the 
eigenvalues problem, with proper boundary conditions assuring continuity at the boundary, see^, with px for Lqq, 
and Vrn(ca) for Lpp. The eigenfunctions of the resonator form a discrete set, which is complete and orthonormal within 
the resonator region. The eigenfunctions of the outside channel region from a continuum indexed by the frequency. 
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w; the index m specifies the asymptotic boundary condition. The functions as the fi\s, form a basis for the 

channel region. 

After the decomposition in the resonators and channel modes, the quantization procedure can be carried out as 
usual: the coefficients multiplying and t'm(w) in the decomposition for A and 11 are promoted to operators with 
proper commutation relations in order to satisfy canonical commutation relations for A and 11. The system-and-bath 
Hamiltonian has, eventually, the form 


n = '^fuvxa{ ctx + '^J duj huj + 


A m 


dui 


WAm(w) a\ Pm{uj) + h.C. 


(14) 


Here, the operators and ax are the bosonic creation and annihilation operators, respectively, associated with the 
resonator mode A; they satisfy = [ 0 ^, 0 :^'] = 0, = Sx,x'- Same for the operators /3m(w) and /?m(w)- 

The matrix Wxmii^) expresses the coupling between the inside and outside functions; explicit expressions can be found 
inM. Note also that in writing Eq. (HI the rotating wave approximations is assumed: we have kept only resonant 
terms (Q;^/3,/3^a) and disregarded nonresonant terms (a^/3^,/3a), which can become important when the widths of 
the modes are comparable to their frequencies. Eventually, since we are interested in the dynamics of the resonator 
modes, we can obtain quantum Langevin equations for the internal modes only: the coupling with the external fields 
gives rise to damping and noise. In the Heisenberg representation, the dynamical equations for the operators ax are: 

ax{t) = -iujxax{t) - tt E [WW^^^,ax'{t)+Fx{t) (15) 

A' 


where Fx{t)^ depending only on the Prn{t) operators, gives the noise term associated to the bath 

Fx{t) = -i f to) (16) 

m 

and we are assuming the Markov approximation: the coefficient Wxm(uj) are supposed to be independent of frequency 
over a sufficiently large frequency band around the atomic frequency^. 

The result of Eq. (1151) generalizes to the case of many modes the Langevin equation obtained for the damping of 
one mode in a cavity, see, for example^. There are two main differences with the equation known in laser theory: the 
noise operators have also non-diagonal correlations, (Fj^Fy) ^ 5xx', and the mode operators ax are not independent 
oscillators, as in standard laser theory, but are coupled through the damping matrix Standard theory is 

recovered in the limit of small leakages, i.e., small coupling with the outside reservoir. 


Langevin equations for mode in random lasers 

Having summarized the main results known for the dumping of multimodes fields, we are now interested in consider¬ 
ing the optical resonator filled with an active medium to actually describe a RL. The medium is usually represented by 
a large number of homogeneously broadened two-level atoms whose density is indicated with p(r). The atomic transi¬ 
tion, uja, is supposed to be much larger than the broadening of the non-diagonal elements of the atomic density-matrix, 
'y±—. In the Heisenberg representation, the atom-field dynamics is described through: 

Aa =-iwAttA - y] 7A/xQf^ + drgl{r)a^{r) + Fx (17) 

(T_(r) = -( 7 _l -biwa)o--(i’) + 2'^gx{r)a:,{r)ax + F_{r) (18) 

A 

^zir) = 7|| {Sp{r) - cr^(r)) - (fflWalcr-W + hx.) -b F^(r) (19) 

A 

where tr^ = \e){g\ and tT_ = \g){e\ are the atomic raising and lowering operator, respectively. |e) and |g) represent 
the excited and ground states of the atom; = |e)(e| — \g){g\j is the inversion operator, 7 || its decay rate. The atom 




density is p(r), while S indicates an external pump injecting atoms in the excited level. In the dipole approximation, 
which we assume to be valid, the atom-field couplings are: 


3 a (r) 


_UJaP_ 

y'2heoU}\ 


MA(r) 


( 20 ) 


being p the atomic dipole matrix element: p oc (erg). To shorten the notation, we have replaced the coupling matrix 
of Eq. (fTKll with 7 ^^ = tt [ITVF’^] 

Semiclassical laser theory consists in neglecting noise terms and replacing operators with their expectation values. 
We consider laser media for which the characteristic times of atomic pump and loss are much shorter than the lifetimes 
of the photons in the resonator. In this case the atomic variables can be removed obtaining nonlinear equations for 
the field alone, equivalent to the master equation for standard lasers, Eq. ©• The standard procedure is to consider 
an expansion in mode amplitudes, ay. one starts neglecting the quadratic term in Eq. (dD, obtaining the zero-order 
approximation, which, replaced in Eq. (1181) gives the first order approximation that replaced back in Eq. m gives the 
second order approximation and so on. We will limit to third order terms; subsequent terms may become important 
far above threshold and, from a statistical mechanical point of view, they are not expected to change the universality 
class of the transition, see, e.g.— . Eventually, the Langevin equations for the mode amplitude take a simple form in 
the basis of the slow amplitude modes. For these modes Ok we make the ansat a^^i^^ : 


a\ Axkakit), with auit) =ak{t)e ( 21 ) 

k 

where the amplitudes Okit) evolve on time scales much larger than the oscillation period These modes play the 
role of the normal modes in a Fabry-Perot resonator for a standard laser, see Eq. and their nontrivial structure 
is a result of the complex multiple-scattering that traps the light in this case. In analogy with the standard case, cf. 
Eq. (|3]), the complex Langevin equations for the amplitude ak{t) can be written then as 

rfH 

«(«) = (22) 

with 

k|FMC(k) ■ k|FMC(k) 


where the sums are restricted to terms that meet the EMC which, for generic 2n interacting modes, reads 


FMC(k) = FMC(ki,..., ka^) : \uJki - H-h “ ^^fc2n I < 7 


(24) 


with 7 being the typical linewidth of the mode. The linewidth can be determined from a full quantum theory, in the 
following it will be supposed as a parameter. The couplings are - S and g^X^Hk^ = ‘^^G^k^k. 2 k 3 ky 

being ^^k^ = and 


P(i-) 9k* (r) 3^ (!■) ’ 


G'iti2fe3fc4 / drpir)g^;ir)g^^{r)gg*{r)g^^{r) (25) 


where where the coefficients proportionality slightly depending on the frequencies involved^, gjt = (A ^)*j, 3 ;i 

and g^ = J2fi^uk9u (cf- Eq- (t20l) 7 Hence, we see the main modifications brought by random lasers: the openness 
of the cavity introduces the linear interaction ^kxk^i due to the coupling of all the modes to the same external bath, 
while the disordered structure of the modes establishes nontrival, positive and negative, interactions gI^\ due to the 
coupling of spatial overlapping modes with the same gain medium. 

Up to now, in the semiclassical approach, we have completely disregarded the effect of noise. However, for a complete 
statistical mechanical description, noise terms must be taken into account providing an effecting bath-temperature 
describing random oscillations. We have seen before that different sources of noise, F\,Fz and F_, occur in Eqs. 
(113-(HU); and, for open cavities, the noise terms, F\, have non-diagonal correlations. The decomposition in the slow 
amplitude modes, Eq. (ED, affects also the noise terms^; this decomposition is by no means unique^ and this 

freedom may be used in order to build a mode basis where the noises do not have non-diagonal correlations. In 

( 2 ') 

general this could be done at the cost to have a non-diagonal linear interactions gk.^k:i ■ study, we 

will consider a situation in which such decomposition is possible and, then, that the noise is white and uncorrelated as 
in the standard laser case, see Eq. (I3|). In this situation, the same theoretical framework adopted for SLD, cf. section 
m can be preserved for general multimode systems: in the limit of small dispersion the evolution becomes Hamiltonian 
and, assuming the spherical constraint £ = Eq, the gain can be considered stationary. This sets up the investigation 
of random lasers through a statistical mechanics study of the Hamiltonian Eq. (1231) . 
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III. A GLASSY LASER 


Among the most singular aspects of RLs is that, at least for some photonic systems composed by a large number 
of modes, a complex behavior in its temporal and spectral response is observed: if there is no specific frequency that 
dominates the others, the narrow emission spikes in the spectra can change frequency from one excitation pulse to 
another with emission spectra that appear different from shot-to-shot’^^n^^^. In these cases the scattering particles 
and all the other external conditions are kept perfectly constant, so these differences are only due to the spontaneous 
emission from which the RL starts at each shot. In these conditions it is observed that the intensity distribution is not 
Gaussian, but rather of the Levy typei^dZ. This is true close to the lasing threshold, whereas far below and far above 
the threshold the statistics remains Gaussian. The shot-to-shot fluctuations provide also a possible explanation for 
the absence of narrow spikes in some experimental studies: if several emission shots are averaged over, with the laser 
being in the above-mentioned complex regime, the spikes can be averaged ould^. Such peculiar behavior demands a 
theoretical explanation. 

A major benefit can be obtained by an innovative point of view based on a statistical mechanics method, analogous 
to the SLD for mode locking lasers described in section H] In this case the theory can profit of the great progress 
made in the last decades in the theory of complex systems. The peculiarity of this approach is that it does not try to 
predict the properties of a given complex system, since “its behavior crucially depends on the details of the system”—, 
but the probability distribution of those properties for systems belonging to a given class. In this case, one is not 
interested in the exact values of the couplings in the Hamiltonian Eq. (I23II in a specific setting but only in their 
general statistical properties. In particular, we point out as in the following the first two moments of the coupling 
probability distributions define alone the thermodynamic behavior of the system. 

In particular, the replica method^ has been adopted to study photonic system based on the Hamiltonian Eq. (1231) . 
In this approach, to probe the multi-state thermodynamic phase space, one considers n copies of the system with 
exactly the same set of disordered couplings (replicas) and evaluates the disorder averaged partition function Z" of 
the replicated system. A continuation to real n is then considered to evaluate the free energy F 


F = -T log Zj = lim -. (26) 

n-i-O n 

As the size of the system grows sufficiently large, the sample-to-sample fluctuations of the free energy die out so that 
it becomes independent of disorder, i.e., it is self- averaging. In the MET one is finally able to get rid of “spatial” 
dependency (corresponding to the frequencies in the photonic case) and express the free energy as a function only of 
the (complex, replica independent) parameter of global coherenee m = analogous to the magnetization 

for spin models, and the (real) generalized overlap matrices 


Qab 


1 

£ 





1 

R,b = ^J2^[a!a'l] . 

^ 1^1 


(27) 


Eor N —>• oo, following the Parisi ansatz, the physical value of the matrices follows from the extremization of the free 
energy functional in the space of the 77.-step RSB matrices^. 

A crucial point is that the mean field theory is exact in realistic physical situations for RLs. This happens if 
the probability distribution of the couplings is the same for all the mode couples (fci, ^ 2 ) and all the mode tetrads 
(fci, k 2 , fca, ki). Such circumstance is, e.g., the one of a RL that presents extended modes with a narrow-bandwidth. 
The last condition assures that \uji — a ;2 + ^3 — a; 4 | < Suj, being Soj the typical linewidth of the lasing modes, so 
that the EMC is always trivially satished. Such physical circumstance is not possible in more traditional applications 
of the theory of complex system, like spin glasses and structural glasses, where the real system is inevitably short 
range. Indeed the nature of the glass phase of short range spin glasses is a debated problem (see, e.g., Ref<^): in 
this case the phenomenology of the system is strongly affected by activated processes, such as nucleation and barrier 
crossing, that are, instead, negligible in the mean field limit. Nevertheless, the degree of localization in the different 
experimental realizations of RLs ranges of several orders of magnitude, from strongly localized to modes extended all 
over the cavity: RLs may then provide an interesting cross-over between the mean field and short range systems that 
might be crucial for the understanding of glass phase in finite dimensions. 

The deep connections between theory of complex systems and disordered photonics has been put forward in the 
works of Angelani et ah ^^i^^ where they first show, using a standard approach for quenched disordered systems based 
on the replica method, that the competition for the available gain of a large number of random modes can lead to a 
behavior similar to that of a glass transition^i^i^: at the lasing transition the photon gas presents an exponential 
number of metastable states (i.e., it has a nonzero complexity) corresponding to different mode locking processes in a 
RL. In this case the leading mechanism for the non-deterministic activation of the modes, underlying the shot-to-shot 
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fluctuations, is identified in frustration and the induced presence of an exponential number of metastable states^S^. All 
the rich phenomenology of glasses, like aging, memory effects and history dependent responses^, is freshly predicted 
for nonlinear photonic systems. In this case the relevant experimental quantities would be, in general, the whole 
correlation functions in time domain 

= {ai{T + t)aiiT))r = + (28) 

where the average is over the time origin r. They note that in molecular glasses the typical microscopic time scales 
are ~ 10“^^ s while for the photonic dynamics are ^ 10“^"^ s (see, e.g., Refi^), so one could earn orders of magnitude 
to study and test theories of glass transition^. However, at our knowledge, the measure of the phase correlations 
required by Eq. (I28p has not been achievable so far, the main complication being the typical low intensity of the RL 
emission, usually not intense enough to successfully use techniques based on second-harmonic generation. 

The analysis has been completed in the subsequent works^L^ obtaining the whole phase diagram for an arbitrary 
degree of disorder, bridging the results of SLD for passive mode locking (cf. section |T| to the aforementioned of a 
completely disordered amplifying medium. In particular, the “glassy” RL transition is present also for a negative value 
of the mean nonlinear interaction. In the thermodynamic construction, negative values correspond to the absence of 
a saturable absorber, cf. Eq. (Cl), and, then, this would imply that mode locking processes in RLs are achieved even 
without a saturable absorber. In Refi^, the same model is also analyzed using the cavity method. We stress that this 
technique could allow to include also the EMC in the description when used on suited graph (cf. section |I]). Work in 
this direction is in progress. 


Results for the general Hamiltonian model 


All the above mentioned results for the statistical approach to RLs are obtained from the same mean-field scheme 
discussed in section|TT]and leading to Eq. (I23L but using two additional approximations: the strong cavity limit and the 
quenched amplitude approximation. The first considers the situation where the leakages from the cavity are very small, 
so that their effects on the linear coupling and on the noise are negligible; they, hence, only contemplate the nonlinear 
term in Eq. The latter assumes a time scale separation between the evolution of “slow” amplitudes and “fast” 

phases of the modes; in this case, they consider an effective Hamiltonian for the phases alone while the amplitudes are 
fixed and have, then, just the role to renormalize the random couplings: = lofcil lofcal kfesi Wkil ■ gkik 2 Hki- 

In section U we have discussed as how the mode locking in standard lasers is achieved as a nontrival locking of the 
mode phases (the so-called “phase waves”, cf. Fig. |I(a) I. In this mechanism the amplitudes play a minor role (de 
facto the same effect is observed also for the “quenched amplitude” version of the model, see the results for the XY 
model in RefJ^). Nevertheless, even in the ordered case, the amplitudes have a nontrival dynamics that leads to 
experimental predictionsiS, cf. also Fig. |l(b)[ The influence of mode amplitude dynamics is even more relevant in 
the case of RLs, where data about phase correlations are hard to obtain (and not available so far, to our knowledge). 
Moreover, the amplitudes are known to have, at least in some experimental configuration, a nontrival behavior, as 
hinted by the shot-to-shot fluctuations of the intensity spectra. 

To include the role of the leakages and the evolution of the mode amplitudes, a series of recent work a^^i^^i^° have 
considered the application of the replica method to the general Hamiltonian Eq. (E51) . This includes the possible 
presence of a linear term and the employment of both the phase and the amplitude of the modes as degrees of freedom 
of the problem. The EMC Eq. (IMl) is, instead, still neglected to obtain a model for which the mean field theory is 
exact. The full inclusion of the EMC requires the use of different techniques, as discussed for the mode locking laser 
in section m 

In the mean-field theory the couplings 3^kxk2k3ki taken as independent identically distributed random 

variables. Without loss of generality, the probability distributions are considered Gaussian with {p = 2,4) 


(Ah,) 


(X exp 




-j, 


(p) 


21 


2cr2 


with 


P'-Jn 


2 _ 

P 2NP-^ 


J, 


(p) 


0 


j(.p) _ 

0 iVP-1 


(29) 


The scaling of ap and with N assures the extensivity of the Hamiltonian. To have a direct interpretation in terms 
of photonics quantities it is also useful to define the photonic parameters 


—aoJo , 

t2_2 t2 

J ^ —OL iJ ^ 


4^^ =(I - ao)Jo , 
Jl =(I - afj^ , 


Rj=^. 


iP=ev^. 


(30) 

(31) 
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a 

(c) Phase diagram in 'P vs q for a high degree of 
disorder Rj (solutions with m = 0) with static (solid) 
and dynamic (dashed) lines. The dynamic line is 
where the dynamics arrests typical of a RFOT occurs 
and an exponential number of metastable states 
appear. On the solid lines between IW, PLW and 
RL(FRSB) the transition is continuous. 



(d) Whole phase diagram in the photonics parameters 
for a = ao- The solid (dashed) red lines correspond to 
continuous (discontinuous) IW-ML transition. The 
blue surface is the RL-ML transition, the orange 
surface the PLW-RL transition and the green surface 
the IW-PLW transition. The two black lines mark the 
intersection between the orange-blue and green-red 
surfaces, respectively. 


The parameters Jq and J fix, respectively, the cumulative strength of the ordered and disordered contribution to the 
Hamiltonian while cto and a the strength of nonlinearity in the ordered and disordered parts. The parameter Rj is 
the degree of disorder and V is the pumping rate, equivalent to the one defined for the mode locking laser in section H] 
The whole phase diagram in terms of the order parameters m, Q and R (cf. Eq. (071) 1 is obtained in Refsj ^^i^^ 
varying the previously mentioned photonics parameters. A crucial point is that the simplest ansatz of assuming 
Qab = Q for all a b, i.e., of assuming that all replicas are equivalent, is not consistent in the whole phase space. 
Therefore, one must allow for a Replica Symmetry Breaking (RSB) to construct the solution, so that the elements 
of Qab and i?ab can take different values and the order parameter becomes their probability distribution. In this 
situation, then, identical copies of the system show different amplitude equilibrium configurations, as the ergodicity 
is broken in many distinct states. 

Four different photonic regimes are found varying parameters V, a, ao and Rj: 

• Incoherent Wave (IW): replica symmetric solution with all order parameters equal to zero. The modes oscillate 
incoherently and the light is emitted in the form of a continuous wave. At low pumping this is the only solution. 
It corresponds to the paramagnetic phase in spin models. 

• Phase Locking Wave (PLW): all order parameters vanish but the diagonal elements of i?ab, so that (ct^) ^ (r^) 
in this phase; this signals a partial locking of the mode phases. This regime occurs in the region of the phase 
space intermediate between the IW and RL regimes. This thermodynamic phase has not equivalent in other 
spin models, since it follows directly from the peculiar nature of the degrees of freedom, that consist of both a 
phase and a magnitude, a combination of the well studied XY and spherical spin models. 

• Mode Locking Laser (ML): solution with m ^ 0, with or without replica symmetry breaking. The modes 
oscillate coherently with the same phase and the light is emitted in form of short pulses. It is the only regime 
at high pumping when the degree of disorder Rj is small enough. It corresponds to the ferromagnetic phase in 
spin models. 

• Random Laser (RL): the modes do not oscillate coherently trivially (i.e., all with the same phase), so that 
m = 0, but the overlap matrices Qab and Rab have a nontrivial structure. It is the only phase in the high 
pumping limit for large degree of disorder Rj. It corresponds to the spin glass phase in spin models. 


The effect of the cavity leakages can be contemplated varying the strength of nonlinearity a and ao. In particular, 
different replica symmetry breaking solutions RL regimes are possible varying a. Consider, e.g., the situation shown 
in Fig. |l(c) of a large degree of disorder Rj, so that Jq -C J and the system is highly frustrated (the scenario is 


stable till the degree of disorder is so low that the ML phase appears). For high a (closed cavity limit), the RL 
regime is IRSB. In this case the transition from PLW to RL is a random first order transition (RFOT) with the 
phenomenology typical of glass transitional: a jump is present in the order parameters Q and R but the internal 
energy remains continuous. For a < ani ~ (3 — I.76382e)(3 — I.03703e^)“^, a continuous part must be considered 
in the distribution of the values for the elements of the overlap matrices. The stable thermodynamic phase is, then, 
first 1-FRSB and eventually a full replica symmetry breaking (FRSB) stated. In realistic optical systems the 2-body 
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Rj Rj 


Figure 1: Phase diagram for e = 1 and a = ao = 1 corresponding to a closed cavity (left) and a = ao = 0.5 corresponding 
to a highly opened cavity (right). The solid (dashed) line between IW and ML corresponds to a continuous (discontinuous) 
transition. The dashed line between PLW and RL in (a) indicates the dynamic transition line with finite complexity. 


interaction is usually not dominant above the lasing threshold, however, there can be systems where the damping due 
to the openness of the cavity is strong enough to compete with the non-linearity and a FRSB state may emerge. In 
this case the transition turns out to be continuous in the order parameters: the overlap is zero at the critical point 
and grows continuously as the power is increased above threshold. 

The whole phase diagram in terms of the four optical regimes is shown in Fig. |l(d)[ and two sections for in Fig. [TJ 
for every degree of disorder and a = ag- Consider the common experimental situation with an increasing pumping 
rate V at fixed Rj, a and ag. Then: 


• for Rj not too large a transition between the IW and the ML regimes is observed increasing the pumping. The 
transition is robust with respect to the introduction of small disorder; 

• for systems with intermediate disorder, the high pump regime remains the ordered ML regime, but the interme¬ 
diate, partially coherent, PLW regime appears between ML and IW; 

• for large i?j, a further transition from ML to RL is observed at high V. Moreover, if Rj exceeds a threshold 
the ML disappears and the only high pumping phase remains the RL. 


In general, the value of ag affects the transition toward the ML regime: for high ag the transition is discontinuous. 
On the contrary, if ag is low, there are regions in the phase diagram where the transition is continuous. The value 
of a controls the transition to the RL regime. For a > ani ( — 0.6297 for e = 1) the tr ansition is toward a RL via a 
RFOT, with a region of finite complexity antecedent the transition (cf. also Fig. 1(c)) where the photonic glass has 
an exponential number of metastable states corresponding to different mode locking processes. 

As we aforementioned, the measure of the phase correlations, required not only by Eq. (l28l) but also for the static 
order parameters Q and R Eq. (E3, is hard because of the low intensity of RL emission and has not been achieved 
so far. Only data for the intensity spectrum are usually available. We note that the lack of a direct experimental 
measure of the overlap probability distribution is common to the original systems for which replica theory was first 
developed, i.e., spin-glasse o^^i^^ , and also to structural glasses, one of the fields of major application of the theory2i“2^. 
The inclusion of the amplitude as degree of freedom in the previous analysis paves the way to a test directly from 
intensity data. Indeed an experimental validation of a random laser-glass connection, and, specifically, of the presence 
of RSB predicted by the theory, has recently been put forward in Refi^S, measuring the overlap between intensity 
fluctuations of the spectral emission. In the model Eq. ( 1 ^ 51 ) . it is possible to define the intensity fluctuations overlap 
matrix 


Cab = 


8A^e2 


N 


k=l 




(32) 


defined in the dominion [0,1]. In the following we consider the symmetrized probability distribution P(|C|), without 
any loss of generality. It can be showit2S that for m = 0 it holds simply Cab = so that the intensity fluctuations 
overlap is nontrivial in the RL regime and the replica symmetry breaking can be studied from intensity spectrum 
data. We note also that instead, for low disorder, at the IW-ML transition Cab does not change at the transition's. 
A typical scenario for the overlap C in an open cavity for high disorder is shown in Eig. [2] at the RL transition 
the distribution becomes nontrivial with more than one value possible, meaning that the replica symmetry is broken. 
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Figure 2: Probability distributions of the C overlap Eq. (1321) for a = oq = 0.4 and Rj = 1.1. Vertical lines represent Dirac 
deltas, whose height is the probability of the argument value. Different regimes are displayed from left to right: the first panel 
correspond to IW (the same distribution is obtained for PLW), the second to RL(l-FRSB) and the last three to RL(IRSB) for 
increasing pumping. 


The behavior shown in Fig. [2] is qualitatively reproduced by the experimental data of Refi^, where the overlap is 
measured between the intensity fluctuations of different shots (assumed as real replicas) in the regime where the 
shot-to-shot fluctuations are observed. These results support, then, the interpretation of the observations as the first 
direct experimental evidence of replica symmetry breaking. 


IV. CONCLUSIONS 

In this paper we have presented the main results of the application of statistical mechanics to nonlinear photonic 
systems, both in standard laser systems with optical cavities and in cavity-less random media. For standard lasers the 
outcome have been quite wide and fruitful, providing a clear and universal description of key phenomena like active 
and passive mode locking and, in particular, the role of the non-linear frequency matching for the onset of pulsed 
emission and the gain narrowing in ultra-fast multimode lasers. 

We have, further, pointed out the crucial role that the statistical mechanic approach can play in the study of 
photonic systems in disordered media, like random lasers. In the latter case, the techniques developed for complex 
systems turn out to be a valid theoretical tool to consider nonperturbatively the concurrent presence of disorder and 
nonlinearity including many aspects not yet established in theoretical optics. In order to have a complete theory for 
light-wave propagation and amplification in disordered active media, e.g., the structure of the lasing modes in an open 
cavity for realistic systems above the lasing threshold has yet to be understood. The presence of a complex structure 
of modes in such systems was revealed already in the first experiments, by the presence of narrow resonances in the 
emission spectrum, emerging as the pumping was increased^S, and supported by later studies of photon statistic a^°’^^ . 
In the course of the last decade, random lasing has been observed in many different kinds of disordered materials, like 
polymer films^S, porous materialsi^, powders^, ceramic a^^’^'^ , clusters^^, colloidal solutions of nanoparticles^, and 
it is now generally regarded as a universal property of optically active disordered structures. However, the physical 
mechanism underlying the constitution of long-living laser modes is, nonetheless, not yet known. Some attempts to 
investigate their formation have been made on the basis of the light localization^^—, which, however, is very difficult 
to accomplish in optics^SiSi. In most materials, indeed, the intensity is spread throughout the sample and the modes 
appear to be extended to the whole system. An alternative explanation is based on the presence of rare long light 
paths in the material!^: the spikes in the emission spectrum are associated with “lucky photons” spontaneously 
emitted in such long paths and acquiring a huge gain, without the need for strong localization. Depending on the 
specific setting, it is expected that both localized and extended modes can, in principle, participate in amplification 
by stimulated emission and which kind of modes is dominant depends, in general, on the delicate balance between 
gain saturation, scattering strength and amplification lengtb^^. 

To discriminate universal and particular characteristics in such composite situation is not an easy task and a 
statistical mechanics approach can help to achieve a classification based on fundamental properties. We have reported 
the main results in this context that, up to now, are limited to mean field theories of effective Hamiltonian models 
of the mode interaction. In particular, we have described the general phase diagrams, ranging from ordered closed 
cavities to disordered open cavities, for the most general model that includes the presence of leakages and the amplitude 
dynamics. An experimental test for replica symmetry breaking at the random lasing threshold can be put forward 
within this theory: the replica symmetry breaking in the intensity fluctuation overlap is shown to be equivalent 
to the one occurring in the standard Parisi overlap (the complex amplitude overlap), allowing for a verification of 
Replica Symmetry Breaking occurrence directly from shot-to-shot fluctuations of emission spectrai^S. The first direct 
experimental evidence, to our knowledge. Our results also hint that a classification can be established between glassy 
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and non-glassy random lasers, depending on the distribution of the disordered nonlinear mode-couplings induced 
by the microscopic structure of wave scatterers and the properties of the optically active material in the candidate 
random laser compound. 

There are several possible developments to the models studied so far for random lasing. The most immediate 
progress would be to go beyond the mean field theory and consider the role of the frequency matching condition 
and the finite spatial extension of the modes. As discussed for the case of standard lasers, this is possible employing 
Monte Carlo simulations or the cavity method on suited graphs. In this respect, we note that the frequency matching 
condition has the same structure present in the Golay-Bernasconi mode l^^i^^ , with discrete Boleean variables, a model 
without quenched disorder but with a complex landscape and a glassy behavior and that has been studied via the 
replica metho d^^’^^ . The most feasible alternative to go beyond mean field is to use numerical techniques like Monte 
Carlo simulations, as accomplished for the ordered case (cf. section |T|. In this way one can also study finite size effects 
that may be relevant in connection with experiments. 

An important counterpart of the analysis here reported is to directly study the dynamics, in its general evolution 
described by the Langevin equations for the electromagnetic modes, without the requirement of equilibrium and 
Hamiltonian dynamics. In this case, in particular, it would be possible to investigate the role of dispersion and of the 
gain dynamics preceding saturation and the onset of the stationary lasing regime. 

A preeminent issue remains the fundamental structure of the lasing modes in cavity-less and disordered media and 
the shape of their effective interaction. A first aspect is that the interactions are not, in general, independent from 
each other, as assumed in the mean field approximation for large number of modes (possible correlations are expected 
to decay with the size of the system). To have access to the values of such interactions in realistic settings is, though, 
not an easy task. We note that, in such an intent, statistical mechanics can again play an important role providing 
the statistical inference techniques to reconstruct the interactions from configuration data (see, e.g. Refi^). 

The benefits obtained so far through statistical physics in the study of laser physics support the use of such 
framework in more setups in future. One peculiar and relevant case might, e.g., be the study of free-electron lasers, 
where the lasing is obtained by high-speed electrons moving in an inhomogeneous magnetic structur o^^i^^ . 

Concluding, we remind that an early motivation for the development of the theory of complex systems was the 
analysis of “certain alloys of ferromagnets and conductors, such as AuFe or CuMg, known as Spin Glasses, [...] 
intriguing and perplexing as no others in the history of solid state physics”^^. The theory actually found quickly 
many areas of application such as structural glasses, neural networks, financial markets, granular materials, the 
immune system, road traffic, flocking in birds or fish, human social communities. After all the above mentioned 
results that have been obtained in the last decades, we feel that also optics can now be entirely added to the list. 
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